A Mechanics Analysis of Carbon Fiber Plain-Woven Thermoset Prepreg during Forming Process Considering Temperature Effect

The preforming quality of carbon fiber plain-woven thermoset prepreg (CFPWTP) is critical to the performance of composite aerospace parts. The deformation ability of the CFPWTP material during preforming is affected by both the fabric woven structure and the resin viscosity, which is different from the dry textile material. Incorrect temperature parameters can enlarge the resin’s viscosity, and high viscosity can inhibit fiber deformation and cause defects. This study proposes an equivalent continuum mechanics model considering its temperature–force behavior. Picture frame tests and axial tensile tests at 15 °C, 30 °C, and 45 °C are conducted to obtain the temperature–stress–strain constitutional equations. By Taylor’s expansion formula and surface fitting method, the constitutive modulus of the material is obtained. Consequently, a saddle-shaped forming simulation is carried out, which is later validated by experiments. Results show that the accuracy of the predicted model is high, with 0.9% of width error and 5.1% of length error separately. Besides, the predicted wrinkles are consistent with the test in fold position and in deformation trend under different temperatures.


Introduction
Advanced composite materials significantly improve performance and reduce aviation parts' weight and cost. Composites, especially the carbon fiber plain-woven thermoset prepreg (CFPWTP), have been increasingly accepted in the aerospace manufacturing industry because of their light weight, high strength, corrosion resistance and formability [1]. When designing composite parts with CFPWTP material, this laminated material usually needs to be preformed layer by layer on the nondevelopable mold surface [2]. At this time, the material will endure great in-plane tensile, compression or shear deformation, which is mainly presented as the interfiber slippage and rotation. Some fiber bundles will be stretched or stacked together if the fiber movement is hindered, and local materials will be bridged or wrinkled [3,4]. Without the correct operation, these defects will seriously damage the aerodynamic shape or interlayer strength of products. Common treatments [5][6][7] include stretching, pressing the defect area, heating or redesigning the cutting size. However, arranging these processes requires firstly determining the deformation characteristics of the CFPWTP material.
The deformation of CFPWTP material is affected by both structure and components. From the perspective of structure, CFPWTP is the product of weaved carbon fiber submerged in resin, which leads to the complex contact interface between these two components. As for the component, CFPWTP inherits the excellent axial mechanical properties and low shearing stiffness [8] of dry woven fabrics as well as the resin's polymer characteristics such as viscoelasticity. Researches [9,10] show that the resin matrix can significantly affect the interply shear resistance of woven fabrics.
In order to obtain material property, many efforts have been devoted through the last decades. Boisse et al. [11,12] and Cao et al. [13] have systematically studied the mechanics of glass fibers with many test methods. They have revealed the characteristic of dry glass fabric materials. Bilisik et al. [14] found that different length-width ratios of samples and the number of pulling fibers influence the shear results significantly. Nosrat-Nezami et al. [15] improved the design of the original picture frame (PF) device and carried out shear tests under different pre-tension conditions. The results show that the shear modulus increases nonlinearly and monotonically with the pre-tension load. Naito et al. [16] experimentally evaluated resin flow mechanism change in prepreg processing during cure including the gelation stage. Wang et al. [17] did work to investigate the interply slipping behaviors to avoid defects such as wrinkling during the manufacturing of large-scale composite parts in hot diaphragm forming.
Besides, some researchers have developed a series of nonlinear constitutive models [18,19] to characterize this material. Peng and Cao [20,21] built the non-orthogonal constitutive model for woven composite. Charmetant et al. [22] established a hyperelastic material model based on continuum mechanics. In addition, based on the woven unit cell, Boisse et al. [23][24][25] created a hyperelastic model which considers the effects of tensile and shear of fiber bundles in mesoscale. Haghi Kashani [26] studied the impact of tension-shear coupling in his model. Chen [27] studied the noncrimp fabrics (NCF) material model and the process optimization method of eliminating wrinkles through a genetic algorithm. Komeili [28] analyzed the fabric model from a multi-scale perspective. In addition, based on the forming theory, three kinds of models have been developed: the dynamic model [29], the discrete model [30] and the continuum model [31]. However, the dynamic and discrete models are not suitable for analyzing the whole formation because of the enormous calculation and lack of material property. Many researchers use the continuous model to simulate the woven fabric composite because it can be embedded in the finite element method (FEM) software [32].
The wrinkle is a typical defect of woven composites [33,34]. Generally, there are two types of wrinkles. 1. In-plane fracture [35]: the material's warp and weft fibers reach the shear locking angle [36][37][38][39], resulting in the inability of the material to continue shear deformation, resulting in the tow buckling and fiber fracture. 2. Out-of-plane bending [40]: the material does not necessarily reach the shear locking angle but folds due to material compression.
There has already been much research on fabrics [41][42][43][44], but few models consider the effect of temperature change on the deformation capacity of materials. For example, it is difficult to move or shear the tow from the woven structure at low temperature due to resin consolidation, while at high temperature when the resin is melted, the above situation is the opposite. Therefore, to characterize the CFPWTP materials, it is necessary to consider both resin's and fiber's properties.
This study establishes an equivalent model of CFPWTP material based on the continuum mechanics for industrial forming processes. The new model is developed from the dry woven-fabric equivalent model but considers the resin's characteristics variety. Then, the picture frame tests and axial tensile tests are carried out at three temperatures to obtain the model's mechanical-temperature characteristic curves. With Taylor's expansion method and fitting of the above results, the constitutive matrix conducive to simulation is obtained. The forming experiments are prepared and compared to examine the simulation's accuracy, and then the wrinkle prediction and the comparison are arranged. The conclusions are prepared in the final section. The structure can be seen in Figure 1.

An Equivalent Continuum Model of CFPWTP Material
This section will introduce the model in three steps. 1. Build a dry fabric model based on the continuum mechanics. 2. Analyze the resin's temperature-viscosity characteristics. 3. Construct the model for characterizing CFPWTP materials.
Before the first step, it is necessary to analyze plain-woven prepreg characteristics for the model building. For example, see Table 1, the parameters of the Toray ® T1100G datasheet of the plain-woven carbon fiber prepreg are offered here. Firstly, through observation, it can be found that the prepreg thickness is very thin. If the composite size is 500 mm × 500 mm, the ratio of length/thickness or width/thickness exceeds 1000. Hence, it is acceptable to set aside the thickness of the model. Because plain weave fabrics 0°and 90° tensile modules are nearing, it is reasonable to consider that the mechanical parameters of the warp and the weft are the same. Researchers treat the plain-weave fabric material as macrolevel periodic materials consisting of a series of mesoscale patterns. This mesoscale pattern is called the representative unit cell (RUC) [46]. For the plain-woven dry fabric, a constitutive model [47] is created as follows. Firstly, acquire the material's characteristic matrix. Secondly, establish the transforming tensor matrix of the global orthogonal coordinate and non-orthogonal coordinate. Thirdly, replace the RUC with the square plate element of FEM, see Figure 2.

An Equivalent Continuum Model of CFPWTP Material
This section will introduce the model in three steps. 1. Build a dry fabric model based on the continuum mechanics. 2. Analyze the resin's temperature-viscosity characteristics. 3. Construct the model for characterizing CFPWTP materials.
Before the first step, it is necessary to analyze plain-woven prepreg characteristics for the model building. For example, see Table 1, the parameters of the Toray ® T1100G datasheet of the plain-woven carbon fiber prepreg are offered here. Firstly, through observation, it can be found that the prepreg thickness is very thin. If the composite size is 500 mm × 500 mm, the ratio of length/thickness or width/thickness exceeds 1000. Hence, it is acceptable to set aside the thickness of the model. Because plain weave fabrics 0 • and 90 • tensile modules are nearing, it is reasonable to consider that the mechanical parameters of the warp and the weft are the same. Researchers treat the plain-weave fabric material as macrolevel periodic materials consisting of a series of mesoscale patterns. This mesoscale pattern is called the representative unit cell (RUC) [46]. For the plain-woven dry fabric, a constitutive model [47] is created as follows. Firstly, acquire the material's characteristic matrix. Secondly, establish the transforming tensor matrix of the global orthogonal coordinate and non-orthogonal coordinate. Thirdly, replace the RUC with the square plate element of FEM, see Figure 2. Because the thickness is negligible, the material element can be regarded as a thin plate. Based on the mechanics of materials and FEM, the thin plate can be treated as a shell element which has the following mechanics matrix ( ) 12  axial  trans  axial11  trans12   21  22  trans  axial  trans21  axial22  33 shear shear can be solved: where the function comes from the PF test and axial tensile test of Figure 3. Because the thickness is negligible, the material element can be regarded as a thin plate. Based on the mechanics of materials and FEM, the thin plate can be treated as a shell element which has the following mechanics matrix [ D(ε)]: where the function comes from the PF test and axial tensile test of Figure 3.
To detect the function's specific expression, it is necessary to establish the relationship between the strain, the stress of the model and the displacement, the force from the tests. The following equations are given in combination with the dimension of the material and device. Specifically, for shear stress δ Shear and the strain ε Shear , Polymers 2022, 14, 2618 5 of 20 where F is the tensile force directly from the test machine, F shear is the solved shear force acting on the material, L mat is the effective length of the cross-shape material, T mat is the thickness of the material, L arm is the arm length of the articulated arm of the PF, and the Path is the displacement of the machine clamp. To detect the function's specific expression, it is necessary to establish the relationship between the strain, the stress of the model and the displacement, the force from the tests. The following equations are given in combination with the dimension of the material and device. Specifically, for shear stress δShear and the strain εShear, where F is the tensile force directly from the test machine, Fshear is the solved shear force acting on the material, Lmat is the effective length of the cross-shape material, Tmat is the thickness of the material, Larm is the arm length of the articulated arm of the PF, and the Path is the displacement of the machine clamp. For tensile stress δAxial and the strain εAxial, they could be solved from the following equations ( ) where F is the tensile force, Wmat is the effective width of the material, Tmat is the thickness of the material, Lmat2 is the effective length of the tensile material, and the Path is the displacement of the machine clamp.

The Orthogonal Coordinate Transforming
Based on the similarity between the square plate element and the RUC pattern, the edge strain vectors of plate elements are used to represent the warp and weft fiber of RUC. Thus, the strain vector of the non-orthogonal coordinate system needs to be transformed into the global orthogonal coordinate for the later FEM analysis. With tensor and continuum mechanics, the relation between strain dε from the orthogonal coordinate and its corresponding dε  from the arbitrary non-orthogonal coordinate is given by For tensile stress δ Axial and the strain ε Axial , they could be solved from the following equations where F is the tensile force, W mat is the effective width of the material, T mat is the thickness of the material, L mat2 is the effective length of the tensile material, and the Path is the displacement of the machine clamp.

The Orthogonal Coordinate Transforming
Based on the similarity between the square plate element and the RUC pattern, the edge strain vectors of plate elements are used to represent the warp and weft fiber of RUC. Thus, the strain vector of the non-orthogonal coordinate system needs to be transformed into the global orthogonal coordinate for the later FEM analysis. With tensor and continuum mechanics, the relation between strain dε from the orthogonal coordinate and its corresponding d ε from the arbitrary non-orthogonal coordinate is given by Using Equation (10) (1) and (10) into Equation (12) can obtain (13) where [D] is the elastic matrix of the orthogonal coordinate system. It is easy to obtain the following relationship

The Analysis of Viscosity Characteristics of the Resin
Different from dry fabrics, the resin viscosity will change while heating, which affects the prepreg's stiffness. In order to investigate the influence of resin, this study used the wp-3011 carbon fiber plain-woven thermoset prepreg (Weihai Guangwei composite material Co., Ltd., Weihai, China) as the test material. The material parameters can be found in Table 2. Three samples are compared to observe the heating effect. Each sample material is a square prepreg of 50 mm × 50 mm. Before the measurement, an incandescent lamp is used to heat the material, and an infrared thermometer is used to measure the sample temperature. Each experiment is conducted three times, and the sample pictures are shown in Figure 4. With the help of the Keyence VHX-600 digital microscope, the observation shows that at 15 • C, the resin is solidified and dispersed (the solidified resin is distributed in the lowest part, almost not connected to each other), and there is less than half the area uncovered by resin. At 30 • C, the semi-cured resin begins to flow, and the resin-uncovered area begins to decrease. At 45 • C, the completely melted resin is evenly distributed on most areas of the fiber.  Because the resin of the prepreg is difficult to separate and measure, this article refers to the resin properties of a similar Toray ® prepreg. From the datasheet of Toray ® G-94M PREPREG SYSTEM, it could be found that before 90 °C, from the data curve of Figure 5, the viscosity and the temperature of the resin show a nearing exponential decline relationship. Because the resin of the prepreg is difficult to separate and measure, this article refers to the resin properties of a similar Toray ® prepreg. From the datasheet of Toray ® G-94M PREPREG SYSTEM, it could be found that before 90 • C, from the data curve of Figure 5, the viscosity and the temperature of the resin show a nearing exponential decline relationship.
Because the resin of the prepreg is difficult to separate and measure, this article refers to the resin properties of a similar Toray ® prepreg. From the datasheet of Toray ® G-94M PREPREG SYSTEM, it could be found that before 90 °C, from the data curve of Figure 5, the viscosity and the temperature of the resin show a nearing exponential decline relationship. From the above observation, it is reasonable to indicate that the temperature affects resin viscosity, which in turn affects yarn sliding ability, and ultimately affects the deformation of CFPWTP. Therefore, the effect of temperature on the viscosity is equivalent to the direct impact of temperature on the CFPWTP deforming ability (stiffness). However, due to the complex interface between fibers, the influence of temperature on CFPWTP seems hard to analyze from the perspective of viscosity alone. In addition, the resin is a polymer material with viscoelastic characteristics. In the same way, temperature also changes its viscoelasticity, which affects the cumulative equivalent elasticity of the CFPWTP. Therefore, the influence of temperature on resin viscoelasticity is also equivalent to the effect of temperature on the CFPWTP deformation ability (stiffness). Moreover, heating affects analysis on CFPWTP in terms of stiffness avoids researching complex contact interfaces. Considering the characteristics, both viscosity and viscoelasticity decrease with heating, and they have similar effects on the deformability of CFPWTP. For simplifying the later analysis, the effect of temperature on the resin is considered entirely on the resin's viscoelasticity.
Based on the above analysis, another assumption is created here: fiber and resin are independent. Therefore, the resin viscoelastic modulus can be theoretically added to the fabric model elastic modulus . Thus, the equivalent Section 2.1 fiber model with the resin characteristics superimposed is considered the model of CFWPTP. The advantage of this method is that the parameter curve of the new CFPWTP model can be obtained through the previous PF and uniaxial tests at specific temperatures. Besides, the From the above observation, it is reasonable to indicate that the temperature affects resin viscosity, which in turn affects yarn sliding ability, and ultimately affects the deformation of CFPWTP. Therefore, the effect of temperature on the viscosity is equivalent to the direct impact of temperature on the CFPWTP deforming ability (stiffness). However, due to the complex interface between fibers, the influence of temperature on CFPWTP seems hard to analyze from the perspective of viscosity alone. In addition, the resin is a polymer material with viscoelastic characteristics. In the same way, temperature also changes its viscoelasticity, which affects the cumulative equivalent elasticity of the CFPWTP. Therefore, the influence of temperature on resin viscoelasticity is also equivalent to the effect of temperature on the CFPWTP deformation ability (stiffness). Moreover, heating affects analysis on CFPWTP in terms of stiffness avoids researching complex contact interfaces. Considering the characteristics, both viscosity and viscoelasticity decrease with heating, and they have similar effects on the deformability of CFPWTP. For simplifying the later analysis, the effect of temperature on the resin is considered entirely on the resin's viscoelasticity.
Based on the above analysis, another assumption is created here: fiber and resin are independent. Therefore, the resin viscoelastic modulus can be theoretically added to the fabric model elastic modulus [ D (ε)]. Thus, the equivalent Section 2.1 fiber model with the resin characteristics superimposed is considered the model of CFWPTP. The advantage of this method is that the parameter curve of the new CFPWTP model can be obtained through the previous PF and uniaxial tests at specific temperatures. Besides, the thermomechanical response envelope surface can then be established with multiple curves at different temperatures, and the fitting equation of the envelope surface is the constitutive equation [ D(ε, T)] of the CFPWTP model, where [ D(ε, T)] is the revised elastic matrix of Equation (14) [ D(ε)] with the temperature parameter T.

Constitutive Model of CFPWTP Material
Because the CFPWTP model is theoretically based on the continuum mechanics. By taking advantage of the differentiable characteristics of the continuum, the [ D(ε, T)] could be expanded with the Taylor's series.
where C and C ij are constants (i, j = positive integer), and the strain ε and the temperature T are the independent variables. Similarly, its elastic matrix can be defined as For calculating the material elasticity matrix D preg (ε, T) of the Equation (16), the following method will be adopted here: (a) Obtain the force-displacement characteristic curve C F-P through the test at a specific temperature, T. (b) Combined with Equations (6)-(9), the characteristic curve C F-P is transformed into the stress-strain curve C δ-ε of the temperature T. (c) Using characteristic curve set {C δ-ε (T1), C δ-ε (T2), C δ-ε (T3), C δ-ε (T4)· · · } to construct the surface S δ-ε-T equation, which characterizes the relationship between temperature T, stress δ and strain ε. (d) Based on Equation (15), calculating the fitting surface equation D ij preg (ε, T) with the nonlinear regression, where i, j = 1, 2, 3.
As for the range of axial tensile curve C F-P , the setting range of the path is from 0 to the distance when the material breaks. Similarly, in the range of the PF tests' curve C F-P , the setting range is from 0 to the distance when the device arms contact. As for the limit of temperature T, the setting value is based on the following considerations: it can be seen from Figure 4 that the resin is approximately solid at 15 • C but melts and evenly distributes at 45 • C. This indicates that the resistance of the resin to the fiber is small enough at 45 • C. In addition, the excessive temperature will aggravate the internal chemical action of the resin and lead to the early curing of the material. Combined with the above analysis, the temperature in this study is set in three values: 15 • C, 30 • C and 45 • C.

Test Device for Material Characteristics
In order to obtain the material property for Equation (15), the PF test and axial tensile test are carried out in this section, and the test devices are shown in Figure 3. The test machine is INSTRON 5982-100 kN. The PF test device arm length L arm is 300 mm, and the size of the cross-shaped material is shown in Figure 6. As for the tensile test, the preparation size of materials is 60 mm × 200 mm. The effective tensile length L mat2 is 100 mm, as shown in Figure 7. The stretching speed of the test machine is always 10 mm/min. In addition, because the fabric material is relatively soft, the material may relax without tension and damage accuracy. The tightening force (20 N) is applied to ensure the material tightening.  In addition, the initial temperature of the test is 15 • C. Due to the large size of the cross-shaped material, controllable hot airflow is used. Its advantage is that the material is surrounded by hot air, which can maintain a uniform temperature field even if the fabric moves.  In addition, the initial temperature of the test is 15 °C. Due to the large size of the cross-shaped material, controllable hot airflow is used. Its advantage is that the material is surrounded by hot air, which can maintain a uniform temperature field even if the fabric moves.
For measuring the heating effect of the PF test, an infrared thermometer is used to detect the temperature values at five points, see Figure 8. The specific values are shown in Table 3. The material temperature error increases as the heating value rises. The error changes from the first 0.1 (15 °C) to the next 0.2 (30 °C) to the last 0.4 (45 °C). It is because the increase in heating leads to easier heat exchange between hot wind and the surrounding air, and it finally leads to an increase in the heat dissipation of materials. In addition, due to the small size of the material in the axial tensile test, only the center area of the test material is measured. The maximum error of the temperature before and after the test at 15 °C, 30 °C and 45 °C are all less than 0.2 °C. Based on the above analysis, the deviation of surface temperature is less than 1 °C. It is supposed that this heating method will not damage test accuracy. For measuring the heating effect of the PF test, an infrared thermometer is used to detect the temperature values at five points, see Figure 8. The specific values are shown in Table 3. The material temperature error increases as the heating value rises. The error changes from the first 0.1 (15 • C) to the next 0.2 (30 • C) to the last 0.4 (45 • C). It is because the increase in heating leads to easier heat exchange between hot wind and the surrounding air, and it finally leads to an increase in the heat dissipation of materials. In addition, due to the small size of the material in the axial tensile test, only the center area of the test material is measured. The maximum error of the temperature before and after the test at 15 • C, 30 • C and 45 • C are all less than 0.2 • C. Based on the above analysis, the deviation of surface temperature is less than 1 • C. It is supposed that this heating method will not damage test accuracy.

Parameters Calculation and Fitting Approaches
In order to obtain the specific form of Equation (15), tensile tests are carried out in

Parameters Calculation and Fitting Approaches
In order to obtain the specific form of Equation (15), tensile tests are carried out in this part to obtain the material characteristic curve. For accuracy and effectiveness, every test is conducted three times, and the average value curves and the error bar are drawn in Figures 9 and 10.

Parameters Calculation and Fitting Approaches
In order to obtain the specific form of Equation (15), tensile tests are carried out in this part to obtain the material characteristic curve. For accuracy and effectiveness, every test is conducted three times, and the average value curves and the error bar are drawn in Figures 9 and 10.  As for the axial tensile test, the curves can be seen in Figure 9. Firstly, the three curves have similar shapes: 1. pre-tensioning stage (0-0.25 mm); 2. straightening stage (0.025-1 mm); 3 fracture stage (1-2 mm). Although the three curves are measured at different temperatures, there is only a small difference in the tensile strength-680 N at 15 °C; 671 N at 30 °C; and 670 N at 45 °C-indicating that the temperature has little effect on the axial tensile strength. The error value is small (<20 N) in the rising stage but becomes large after a fracture (>100 N). Thus, the curve-fitting process is only applied to the curves before fracture. Through comparison, the max amplitude of the error bar of the three curves is about 20 N, due to the fact that the max error value/maximum value (20 N/670 N) of the curve is less than 3%. In this paper, it is believed that the curve can characterize the axial properties of the material As for the shearing test, three shear curves show similarity. From Figure 10, on the As for the axial tensile test, the curves can be seen in Figure 9. Firstly, the three curves have similar shapes: 1. pre-tensioning stage (0-0.25 mm); 2. straightening stage (0.025-1 mm); 3 fracture stage (1-2 mm). Although the three curves are measured at different temperatures, there is only a small difference in the tensile strength-680 N at 15 • C; 671 N at 30 • C; and 670 N at 45 • C-indicating that the temperature has little effect on the axial tensile strength. The error value is small (<20 N) in the rising stage but becomes large after a fracture (>100 N). Thus, the curve-fitting process is only applied to the curves before fracture. Through comparison, the max amplitude of the error bar of the three curves is about 20 N, due to the fact that the max error value/maximum value (20 N/670 N) of the curve is less than 3%. In this paper, it is believed that the curve can characterize the axial properties of the material As for the shearing test, three shear curves show similarity. From Figure 10, on the whole, the three curves are nonlinear and similar to exponential curves. However, with careful analysis, it can be found that the curve contains approximately linear parts (0-100) and fast-rising parts (100-160). From the view of the curves' maximum value at the linear stage (0-100 mm), the max force of 15 • C is approximately 60 N. At the same stage, the 30 • C ultimate value is about 30 N, while the 45 • C value is 20 N. When the material enters the nonlinear stage (100-160), the maximum value of the 15 • C is about 1600 N, which ranks first, and the value at 30 • C is second at about 790 N. The ultimate value of the 45 • C curves is about 300 N at the third level. Because the maximum curve value (300 N) at 45 • C is 18.7% of the curve value (1600 N) at 15 • C, temperature influence on shear deformation is not negligible during analysis. The above material curve characteristics can be explained from the following aspects. In the first nearing linear parts, the prepregs experienced pure shear condition where the warp and weft tows rotate then compact, resulting in the force increasing gradually. In the nonlinear parts (100-160), the fiber rotation is restricted due to shear locking, and the force along the fiber induces a sharp increase. As for the error, in the rising stage, the max error of the 15 • C curves is 30 N, the max value is 31 N for the 30 • C, and the error amplitude is 20 N for the 45 • C. Because the max error value/maximum (31/790) value of the curve is no more than 7%, and the error is evenly distributed on both sides of the curve, these curves are considered to be sufficient to characterize the shear properties of the material.
Because the property curves are available, a polynomial fitting operation is performed in this section, and the R 2 method is used to judge the fitting accuracy. The criterion is that R 2 closer to 1 means better quality. The calculating equations are below. (17) where RSS is the error sum of squares, TSS is the total sum of squares, where n is the number of original data, y i is the i-th data, f i is the i-th predicted value, and y is the mean value of total n data. Table 4 shows that the R 2 values of all the fitted curves are very close to 1, which indicates that the fitted curves can effectively characterize the original curves. Next, to obtain D preg (ε, T) , the following operations are arranged here: 1. Convert the polynomial curve equations at each temperature into a stress-strain relation curve according to the Equations (6)- (9). The curves of the modulus are presented in Figures 11 and 12 , the following operations are arranged here: 1. Convert the polynomial curve equations at each temperature into a stress-strain relation curve according to the Equations (6)- (9). The curves of the modulus are presented in Figures 11 and 12.   Figure 11. Axial stress-strain curves at three temperatures.   , the following operations are arranged here: 1. Convert the polynomial curve equations at each temperature into a stress-strain relation curve according to the Equations (6)- (9). The curves of the modulus are presented in Figures 11 and 12.   Figure 11. Axial stress-strain curves at three temperatures.  2. Arrange the stress-strain fitting curves in the coordinate system in order of the temperature from low to high. Create an envelope surface for the above three curves. This surface is the stress-strain-temperature characteristic surface of CFPWTP.
3. According to Equation (16), the fitting operation is performed on the envelope surface. The shear and axial characteristic envelope surfaces are shown in Figures 13 and 14.
From values are very close to 1, it also indicates that the surface equations have fully included the characteristics of the modulus curves. Besides, the traditional method only tests at one temperature and obtains the fixed elastic modulus of Equation (18). The characteristic modulus obtained in this study is able to change with temperature, which can be used to determine the optimal temperature with broader adaptability.   the characteristics of the modulus curves. Besides, the traditional method only tests at one temperature and obtains the fixed elastic modulus of Equation (18). The characteristic modulus obtained in this study is able to change with temperature, which can be used to determine the optimal temperature with broader adaptability.

Introduction of Simulation and Test
Because the CFPWTP characteristic matrix D preg (ε, T) is available, this section presents a finite element-forming simulation on the Abaqus CAE 2016 -SIMULA TM by Dassault Systèmes ® , Paris, France. The simulation model includes four parts: the punch, the cover plate, the CFPWTP material, and the die mold, see Figure 15. Specifically, the model's element type is S4R and its thickness is 0.2 mm. According to Equation (18), the material properties are written in the VFABRIC subroutines. Besides, the punch is a saddle-shaped mold with double curvature, and the geometric feature of the punch is symmetrical about the X-axis and Y-axis. The top surface is flat and the height of punch is 40 mm, and its side wall is part of a large round barrel (R 2 = 400 mm). The side wall and bottom surface are transited by round corners (R 4 = 10 mm). The front and rear end faces of the punch are an arc with a radius of 40 mm. In particular, the distance between the die and the punch is 1.8 mm. The CFPWTP material's size is 500 mm × 500 mm × 0.2 mm, and the friction coefficient between every two surfaces is set to 0.1. The operation process of simulation is as follows. the die and the punch is 1.8 mm. The CFPWTP material's size is 500 mm × 500 mm × 0.2 mm, and the friction coefficient between every two surfaces is set to 0.1. The operation process of simulation is as follows. 1. Place the material on the die in advance. 2. Place the cover plate on the material, where the angle between the material and the X-axis is 45°. 3. Set the temperature at 15 °C, 30 °C, and 45 °C separately. 4. Apply concentrated forces 8 × 50 N on the cover plate and 2000 N on the punch. The simulation is run on ThinkPad P15v i7-2.60 GHz, and the results are given in Figures 15 and 16.

Forming Simulation and Test
From Figure 15, the deformed material shows an approximated diamond shape from the top, and the material surface without cover plate pressing is full of wrinkles. Because the whole material is symmetrical, only the third quadrant section will be analyzed and discussed. the die and the punch is 1.8 mm. The CFPWTP material's size is 500 mm × 500 mm × 0.2 mm, and the friction coefficient between every two surfaces is set to 0.1. The operation process of simulation is as follows. 1. Place the material on the die in advance. 2. Place the cover plate on the material, where the angle between the material and the X-axis is 45°. 3. Set the temperature at 15 °C, 30 °C, and 45 °C separately. 4. Apply concentrated forces 8 × 50 N on the cover plate and 2000 N on the punch. The simulation is run on ThinkPad P15v i7-2.60 GHz, and the results are given in Figures 15 and 16.

Forming Simulation and Test
From Figure 15, the deformed material shows an approximated diamond shape from the top, and the material surface without cover plate pressing is full of wrinkles. Because the whole material is symmetrical, only the third quadrant section will be analyzed and discussed.

Forming Simulation and Test
From Figure 15, the deformed material shows an approximated diamond shape from the top, and the material surface without cover plate pressing is full of wrinkles. Because the whole material is symmetrical, only the third quadrant section will be analyzed and discussed.
For quantitative analysis, the length and width of deformed material at three temperatures are presented in Table 5. A comparison of these values shows that in the width direction (X direction), the length of the material increases first and then decreases, but the change is very small (<4 mm), basically maintained at 340 mm at three temperatures. However, in the length direction (Y direction), the size of the material decreases when the temperature rises. Specifically, the initial value is almost unchanged at 353 mm, but it drops sharply to 336.2 mm at 45 • C. In addition, there is an uncovered triangular area at the lower-left corner of the mold after stamping. As the temperature increases, the material's deformation intensifies, resulting in the uncovered blue area becoming larger. These phenomena are because the deformation difficulty of the CFPWTP decreases with the heating, resulting in easier deformation of the material. To further analyze the character of the material, the in-plane shear stress SFABRIC12 of the material is analyzed, as shown in Figure 16. From the stress ranges, it can be seen that when T = 15 • C, the amplitude range of SFABRIC12 is the largest (−5.442 × 10 6 , 1.116 × 10 7 ), and the overall distribution is uneven. At 45 • C, the amplitude range of SFABRIC12 is the smallest, the overall distribution is uniform, and the average value is about (−1.667 × 10 5 , 1.667 × 10 5 ). The value range at 30 • C is (−6.97 × 10 5 , 5.37 × 10 5 ) which is between them, and the stress distribution is more uniform than at 15 • C. From the point of view of stress distribution, the stress of the whole surface presents a uniform state at 45 • C. However, at 30 • C, the maximum positive shear stress and the maximum negative shear stress are mainly near the position where the cover plate is pressed. Large stress is much more concentrated at the lower left side of the die mold. At 15 • C, the whole surface shear stress changes sharply, and the shear stress range enlarges significantly. This shows that as the temperature decreases, folds will first appear near the pressed area of the cover plate. As the temperature decreases, the material deformation becomes more difficult, and the folds will further spread to the whole surface. From the above analysis, it can be concluded that the shear stress of the material decreases and tends to be uniform with the temperature rise.

Forming Test and Comparison
To examine the simulation accuracy, the forming test with the 500 mm × 500 mm wp-3011 prepreg is arranged in this part. The punch, cover plate and die mold are made of aluminum alloy. The punch pressing force is 2000 N, and the cover plate is applied with four air cylinders (4 × 100 N). Besides, the die mold is equipped with a temperature controller for preheating purposes. The resin's melting could affect the friction coefficient between the material and the punch and the die surface. For preventing the unpleasant friction coefficient changes, the whole surfaces of the punch, the die and the cover plate are pasted with the smooth Teflon coating to stabilize the friction. Then the test results are as follows.
By drawing a grid on the photos, the size of the material is measured, as shown in Figure 17. The size errors of the simulation and test are calculated, and the data are shown in Table 5 and Figure 18. It can be seen that the material length in the X -direction will increase first and then decrease while heating, but the magnitude (about 342 mm) of the change is small (<4 mm). As for the material length in the Y-direction, the length first decreases slightly from 346.4 mm at 15 • C to 343.4 mm at 30 • C and then decreases substantially from 30 • C to 353.5 mm at 45 • C. In addition, overall, the X err is generally smaller than Y err , in which the maximum error in the X direction is 0.9%, while the maximum error in the Y direction is 5.1%. The possible reason for simulation and test difference is that the mold shape makes the material easier to shrink in the Y direction, so the deformation in the Y direction is more significant, bringing more error. As for the shape, it can be seen from Figure 15 that as the temperature rises, the uncovered area at the lower-left corner of the die becomes larger, which is consistent with the trend predicted by the simulations. Finally, by comparing the test results and simulation, it can be seen from Figure 18 that although there are shape deviations, the maximum error is not more than 5.1%, and the simulation and test results show the same deforming trend. Therefore, the model established in this paper can still meet the needs of following defect prediction.
shape, it can be seen from Figure 15 that as the temperature rises, the uncovered area at the lower-left corner of the die becomes larger, which is consistent with the trend predicted by the simulations. Finally, by comparing the test results and simulation, it can be seen from Figure 18 that although there are shape deviations, the maximum error is not more than 5.1%, and the simulation and test results show the same deforming trend. Therefore, the model established in this paper can still meet the needs of following defect prediction.

Wrinkles Prediction and Comparison
The typical defect of woven composites is the wrinkle. In order to predict the possible defect location with the established model, two criteria of the wrinkles are established. 1. According to the reference [13], the maximum shear angle should not exceed the shear locking angle of 45° (in-plain wrinkle). 2. The normal distance from the mesh node to the die surface should be equal to the material thickness (out-of-plain wrinkles).
Through quantitative analysis, Figure 19 shows that at 15 °C, the shear angle of every element is below the shear lock angle, but there are obvious folds whose max shape, it can be seen from Figure 15 that as the temperature rises, the uncovered area at the lower-left corner of the die becomes larger, which is consistent with the trend predicted by the simulations. Finally, by comparing the test results and simulation, it can be seen from Figure 18 that although there are shape deviations, the maximum error is not more than 5.1%, and the simulation and test results show the same deforming trend. Therefore, the model established in this paper can still meet the needs of following defect prediction.

Wrinkles Prediction and Comparison
The typical defect of woven composites is the wrinkle. In order to predict the possible defect location with the established model, two criteria of the wrinkles are established. 1. According to the reference [13], the maximum shear angle should not exceed the shear locking angle of 45° (in-plain wrinkle). 2. The normal distance from the mesh node to the die surface should be equal to the material thickness (out-of-plain wrinkles).
Through quantitative analysis, Figure 19 shows that at 15 °C, the shear angle of every element is below the shear lock angle, but there are obvious folds whose max

Wrinkles Prediction and Comparison
The typical defect of woven composites is the wrinkle. In order to predict the possible defect location with the established model, two criteria of the wrinkles are established. 1. According to the reference [13], the maximum shear angle should not exceed the shear locking angle of 45 • (in-plain wrinkle). 2. The normal distance from the mesh node to the die surface should be equal to the material thickness (out-of-plain wrinkles).
Through quantitative analysis, Figure 19 shows that at 15 • C, the shear angle of every element is below the shear lock angle, but there are obvious folds whose max normal height is about 2 mm. This type of fold belongs to out-of-plane bending, see red arrows. At 30 • C, the material has in-plane and out-of-plane folds. Among them, the most serious out-of-plane fold height is 0.8 mm, and the maximum shear angle of the in-plane fold is 51 • . At 45 • C, the element fluctuation is lower than that of the 30 • C conditions, and the amplitude is less than 0.4 mm. The largest shear angle is 29 • . Moreover, as for wrinkles distribution, at 15 • C, folds are mainly distributed at the corner of the front end of the die and the transition area of the side wall. At 30 • C, the out of plane folds are distributed in the transition region, and the in-plane folds are distributed at the left top edge of the mold. At 45 • C, the whole surface is smooth, but there are only small folds at the sharp corners of the mold. Compared with the shear stress distribution of SFABRIC12, the wrinkles' location shows consistency with the large shear stress position on the mold.
The above phenomenon happens because of a few reasons. 1. The geometric constraint: the press of the cover plate makes it difficult for the material at this position to bulge substantially, so the folds at this position are generally in-plane folds. In the transition region of the sidewall, the material is in a free state, and the material will turn 90 • in this region, which leads to the region being in a bending state all the time, and it is very easy to produce out of plane folds. 2. The temperature effect: the temperature affects the shear deformation stiffness of materials. When the temperature is high enough, e.g., 45 • C, the shear stiffness and shear stress will decrease, both of which are conducive to the smooth and large deformation of materials without folds. When the temperature is not high enough, e.g., 30 • C, the material may have both in-plane and out-of-plane folds, depending on the situation and constraint of the material. At low temperature, e.g., 15 • C, the material is difficult to deform. At this time, the internal stress of the material needs to be released by the deformation which is prone to causing the out-of-plane defect.
normal height is about 2 mm. This type of fold belongs to out-of-plane bending, see red arrows. At 30 °C, the material has in-plane and out-of-plane folds. Among them, the most serious out-of-plane fold height is 0.8 mm, and the maximum shear angle of the in-plane fold is 51 °. At 45 °C, the element fluctuation is lower than that of the 30 °C conditions, and the amplitude is less than 0.4 mm. The largest shear angle is 29°. Moreover, as for wrinkles distribution, at 15 °C, folds are mainly distributed at the corner of the front end of the die and the transition area of the side wall. At 30 °C, the out of plane folds are distributed in the transition region, and the in-plane folds are distributed at the left top edge of the mold. At 45 °C, the whole surface is smooth, but there are only small folds at the sharp corners of the mold. Compared with the shear stress distribution of SFABRIC12, the wrinkles' location shows consistency with the large shear stress position on the mold. The above phenomenon happens because of a few reasons. 1. The geometric constraint: the press of the cover plate makes it difficult for the material at this position to bulge substantially, so the folds at this position are generally in-plane folds. In the transition region of the sidewall, the material is in a free state, and the material will turn 90° in this region, which leads to the region being in a bending state all the time, and it is very easy to produce out of plane folds. 2. The temperature effect: the temperature affects the shear deformation stiffness of materials. When the temperature is high enough, e.g., 45 °C, the shear stiffness and shear stress will decrease, both of which are conducive to the smooth and large deformation of materials without folds. When the temperature is From the test results, Figure 20 shows that at 15 • C, there are two big wrinkles on the side wall with a height of 5 mm (see red arrows). According to Figure 15, there is also some waviness at the outer edge of the mold with a height of 1 mm. At 30 • C, the side wall folds are reduced to one, and the height is reduced to 3 mm. The waviness in the outer edge of the mold area is aggravated, but the height is less than 1 mm, as shown in Figure 15. At 45 • C, the original folds shrink, and a new small wrinkle in the corner appears, and the height of it is about 1 mm. The waviness in the top cover plate area also disappears. The big folds are mainly distributed in the outer space of the mold.
By comparing the simulation and test results, the simulation is consistent with the test in fold position and deformation trend under different temperatures. It is considered that the prediction effectiveness has been verified. However, it can be found that there is a certain deviation between the model and the test, such as the model failing to predict the waviness at the outer edge of the mold at 15 • C. It is because firstly, in order to simplify the model, only three temperatures are considered in the model; secondly, the coupling effect inside the material is ignored in the model. To decrease the above deviation, the model will add more parameters in future work. some waviness at the outer edge of the mold with a height of 1 mm. At 30 °C, the side wall folds are reduced to one, and the height is reduced to 3 mm. The waviness in the outer edge of the mold area is aggravated, but the height is less than 1 mm, as shown in Figure 15. At 45 °C, the original folds shrink, and a new small wrinkle in the corner appears, and the height of it is about 1 mm. The waviness in the top cover plate area also disappears. The big folds are mainly distributed in the outer space of the mold. By comparing the simulation and test results, the simulation is consistent with the test in fold position and deformation trend under different temperatures. It is considered that the prediction effectiveness has been verified. However, it can be found that there is a certain deviation between the model and the test, such as the model failing to predict the waviness at the outer edge of the mold at 15 °C. It is because firstly, in order to simplify the model, only three temperatures are considered in the model; secondly, the coupling effect inside the material is ignored in the model. To decrease the above deviation, the model will add more parameters in future work.

Conclusions
The fabric material deforming ability is crucial for composite aerospace part manufacturing. An equivalent continuum model is presented for predicting the mechanics of CFPWTP material. Based on the analysis of materials, this paper arranges mechanics tests at 15 °C, 30 °C and 45 °C. By Taylor's expansion formula and the surface fitting method, the surface constitutive equation of the material is obtained. Besides, a forming test is carried out to examine the model's effectiveness. Results prove that 1., the increase in temperature can reduce the deformation stiffness of the material and effectively inhibit the folding of the CFPWTP; 2., the CFPWTP model in this paper contains the temperature characteristic of the resin and has high accuracy-the width error is less than 0.9%, and the length error is less than 5.1%; 3., the fold position predicted by simulation is basically similar to the test results. The models and methods in this paper can provide effective tools for designers. In addition, to make further improvement to the model accuracy, the model will introduce more parameters such as parameter coupling in the follow-up study.

Conclusions
The fabric material deforming ability is crucial for composite aerospace part manufacturing. An equivalent continuum model is presented for predicting the mechanics of CFPWTP material. Based on the analysis of materials, this paper arranges mechanics tests at 15 • C, 30 • C and 45 • C. By Taylor's expansion formula and the surface fitting method, the surface constitutive equation of the material is obtained. Besides, a forming test is carried out to examine the model's effectiveness. Results prove that 1., the increase in temperature can reduce the deformation stiffness of the material and effectively inhibit the folding of the CFPWTP; 2., the CFPWTP model in this paper contains the temperature characteristic of the resin and has high accuracy-the width error is less than 0.9%, and the length error is less than 5.1%; 3., the fold position predicted by simulation is basically similar to the test results. The models and methods in this paper can provide effective tools for designers. In addition, to make further improvement to the model accuracy, the model will introduce more parameters such as parameter coupling in the follow-up study.